Importance of electronic correlations for structural and magnetic properties of the 

iron pnictide superconductor LaFeAsO 
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We present calculations of structural and magnetic properties of the iron-pnictide superconduc- 
tor LaFeAsO including electron-electron correlations. For this purpose we apply a fully charge 
self-consistent combination of Density-Functional Theory with the Dynamical Mean-Field theory, 
allowing for the calculation of total energies. We find that the inclusion of correlation effects gives 
a good agreement of the Arsenic z position with experimental data even in the paramagnetic (high- 
temperature) phase. Going to low temperatures, we study the formation of the ordered moment in 
the striped spin-density- wave phase, yielding an ordered moment of about 0.60 /xb, again in good 
agreement with experiments. This shows that the inclusion of correlation effects improves both 
structural and magnetic properties of LaFeAsO at the same time. 

PACS numbers: 71. 15. Mb, 71.20.Be, 74.70.-b 
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I. INTRODUCTION 

Since the discovery of high-temperature superconduc- 
tivity in iron-based compounds^ a lot of research has been 
dedicated to this fascinating class of materials. On the 
theoretical side, calculations within density-functional 
theory (DFT), often performed within the local den- 
sity approximation (LDA), could reproduce a variety 
of properties, but failed in the quantitative description 
of other features, such as the mass renormalization of 
the predominately iron-d quasi-particlcs, which could be 
improved by the inclusion of correlation effects for the 
Fe 3d electrons^— A very puzzling mystery shows up 
for the compounds exhibiting long-range spin-density- 
wave (SDW) magnetic ordering at low temperatures. Al- 
though the spin-pattern and ordering vectors were well 
predicted by DFT, it was soon realized that there is a 
big discrepancy between the magnitude of the measured 
magnetic moments with theoretical predictions of spin- 
polarized DFT calculations. For instance for LaFeAsO, 
early experimental data pointed to a very small ordered 
moment in the range of range of around 0.3 ^bt^ al- 
though recent measurements indicated a somewhat larger 
moment of 0.63 in LaFeAsO^ On the other hand, 
DFT calculations using the experimental crystal struc- 
ture always gives large values between 1.7 and more than 
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There is also a strong connection between the value of 
the ordered moment and details of the crystal structure. 
As stated above, the ordered moment turns out to be 
too large in DFT calculations, but in this magnetic case, 
structural optimization of the ^-position of the arsenic 
ions reproduces well the experimental position^ On 
the other hand, non-magnetic DFT calculations, which 
should correspond to the paramagnetic high-temperature 
phase, gives a too short Fe-As distance with drastic in- 
fluence on the low-energy electronic structured Exper- 



imentally, the As z position hardly changes across the 
magnetic transitions^ a fact that is hard to reconcile 
within DFT calculations, since the optimized internal 
structural parameters differ significantly between mag- 
netic and non-magnetic calculations. The correct de- 
scription of the equilibrium structure is particularly im- 
portant for cases where the forces on the ions are impor- 
tant, e.g. phonon calculations. 

There were several attempts to improve over simple 
DFT calculations. Concerning the ordered moment, 
Yildirim et alr^ performed fixed moment DFT calcula- 
tions in order to study the stability of magnetic ordering 
patterns. Attempts to include correlation effects by per- 
forming LDA+U or GGA+U calculations were not suc- 
cessful. It has been shownii that the magnetic moment 
even increases with increasing U, and even for small in- 
teraction values of U f» 1 eV the topology of the Fermi 
surfaces is changed drastically, incompatible with exper- 
iments. However, a reduction of the magnetic moment 
could be found in LDA+U calculation using an effective 
negative interaction parameter Particularly promising 
are approaches using many-body techniques to include 
electronic correlation effects. Using a combination of 
DFT with the dynamical mean-field theory (DMFT) a 
significant reduction of the magnetic moment could be 
found for BaFe2As2 j 1 ^ variational Monte Carlo gave sim- 
ilar results also for other materials,— consistent with a 
recent comprehensive LDA+DMFT study for a variety 
of pnictide and chalcogcnidc materials i2i A general ar- 
gument is that quantum fluctuations hinder a large in- 
stantaneous iron moment from orderingi 1 ^— 

Regarding the combination of structural and mag- 
netic properties, one proposal for a better description 
is to combine in a sophisticated way magnetic and non- 
magnetic DFT calculations.— This approach has been 
used to study the electron-phonon interaction in iron- 
pnictidc superconductors. Improved structural optimiza- 
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tion has been performed using a combination of DFT 
with Gutzwiller wave function techniques^ where the 
values of the interaction parameters where fitted to give 
the correct As height above the Fe plane. 

The motivation for this article is to show that the in- 
clusion of correlation effects by LDA+DMFT for the de- 
scription of LaFeAsO improves substantially the agree- 
ment of both the As z position as well as the ordered mag- 
netic moment between theory and experiment within one 
set of ab-initio calculated interaction parameters, which 
are determined with the constrained Random Phase Ap- 
proximation (cRPA) i 27 i 28 This kind of investigation has 
not been done so far, since theoretical studies including 
strong electron-correlations were focused on the calcula- 
tion of either magneti o 19 i 21 i 22 or structural properties^ 

A consistent approach to total energy calculations and 
structural optimization within LDA+DMFT (the As z 
position in LaFeAsO in the present work) requires self- 
consistency over the charge density. LDA+DMFT is 
often employed within the so-called "one-shot" scheme, 
where the one-electron part of the Hamiltonian obtained 
from the band-structure LDA part is not updated dur- 
ing subsequent DMFT calculations. However, correlation 
effects will in general induce a certain redistribution of 
the charge density, which in turn leads to different Kohn- 
Sham potential and one-particle part of the Hamiltonian. 
The correlation-induced changes in the charge density 
and one-electron potential will also affect the electron- 
nuclei, Hartree, and exchange-correlation contributions 
to the LDA+DMFT total energy. Moreover, in some 
systems (e.g. cerium oxides^ 9 -) the charge density self- 
consistency has been demonstrated to be important for 
spectral properties as well. 

The paper is organized as follows. In Sect. [H] we intro- 
duce the full charge-self consistent LDA+DMFT method, 
followed by Sect. IIIII where we present results for the 
LaFeAsO system. We draw our conclusions in Sect. IIV1 
which are followed by App. [X] with a more detailed dis- 
cussion of the influence of the full charge self-consistency 
on the single-particle spectra. 



II. METHODS 

For the present study, we use a further development of 
a previously introduced LDA+DMFT implementation, 
Ref. H, which is based on the full potential (linearized) 
augmented plane wave method as implemented in the 
Wicn2K packaged Our task of optimizing the arsenic 
ion position necessitates rather accurate calculations of 
the total energy, which, as explained in the Sec. |H re- 
quire a LDA+DMFT scheme fully self-consistent in the 
charge density. The implementation of full charge self- 
consistency is currently a topic of high interest, and sev- 
eral schemes have been implemented recently . 29 ' 31 " — 

Within the projective technique for formation of the 
correlated orbitals, Ref. [f| we use the Kohn-Sham 
(KS) states within a chosen energy window W to form 



Wannier-like functions that are treated as correlated or- 
bitals. In the present work, we use an energy window 
from — 6.8eV to 2.8 eV, spanning the range of Fe-d as 
well as As-p and O-p states, giving a total number of 22 
bands inside the window. Onsite interactions were then 
applied to the five Fe-d orbitals. The very same projec- 
tion scheme has already been used in Ref. [5|. Solving the 
corresponding single-site quantum impurity problem pro- 
duces the local self-energy within the correlated orbitals 
basis set, which is then upfolded into the lattice self- 
energy (k, iuj n ), where lu„ are Matsubara fermionic 
frequencies . The lattice self-energy S(k, iuj n ) is gener- 
ally non-diagonal in the subspace of the KS eigenstates 
{+} (y G W) leading to a non-diagonal lattice Green's 
function within W and to the corresponding density ma- 
trix: 

Nt> =^G^(k^„)e i; "" 0+ (1) 

n 

being also non-diagonal. The charge density distribution 
in the real space is then calculated from the density ma- 
trix N* v , as follows: 

Pdmft (r) = p ow (r) + (r\^ku)N^,(^ kv/ \r), (2) 

where ^ku are the KS eigenstates within the energy 
window W, p ow {v) is the contribution from states out- 
side W. By substituting into @ the expansion of the 
KS eigenstates within the linear augmented plain-wave 
(LAPW) basis set one derives formulas for the charge 
density within the muffin-tin (MT) spheres and in the 
interstitial. These formulate are generalizations of the 
standard LAPW expressions to the case of a density ma- 
trix non-diagonal in the space of KS states. As in the 
standard case, the charge density within the MT spheres 
is expressed through radial solutions (and their energy 
derivatives) of the corresponding Schrodingcr equation. 
In the interstitial it is expressed through plain waves. 
The derivation and relevant formulas for each case are 
given in Appendix |B| 

The LDA+DMFT total energy reads 3 ^ 

E =Eun + E c [pdmft] + Eh[pdmft]+ 
E xc [pdmft] + (Hu) — Ebc, 

where the corresponding contributions in the right- 
hand side are the kinetic, crystal (electron-nuclei and 
nuclei- nuclei) , Hartree, exchange-correlation, Hubbard 
and double-counting correction terms, respectively. The 
second, third and four terms arc evaluated in accordance 
to the standard DFT-LDA expressions but with the up- 
dated LDA+DMFT charge density ©. The kinetic en- 
ergy contribution reads 

Ekin = E band - / dYV K s(Y)p DMFT {Y), (4) 
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where the Kohn-Sham potential vks corresponds to the 
LDA+DMFT charge density pdmft an d the band en- 
ergy contribution E^and is 



III. RESULTS 
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where H^ s is the one-particle (Kohn-Sham) part of the 
Hamiltonian, are its eigenstates with v £ W, and 
Eband i s the sum over the occupied Kohn-Sham eigen- 
states laying outside of the window W. 

Finally, the Hubbard term (Hjj) was evaluated 
in accordance with the Migdal formula (Hjj) = 
■iTr (S(iw)G(iaj)), where and G{iuj)) are the im- 

purity self-energy and Green function, respectively. 

For the solution of the quantum impurity problem we 
apply the continuous-time quantum Monte Carlo method 
in the strong-coupling formulation^ Restricting our- 
selves to density-density interactions only, as in Ref. H, 
we are able to perform calculations down to temperatures 
as low as T = 77 K with reasonable numerical effort and 
without further approximations. For the calculation of 
total energies, high-quality numerical data is necessary. 
In order to get an estimate of the statistical error on the 
total energy, we perform several further iterations (order 
10) after self-consistency is reached, yielding an estimate 
for the standard deviation. 

Interaction parameters have been calculated previously 
within cRPA^, yielding an average Coulomb interaction 
of U = 2.7 eV and Hund's exchange of J = 0.8 eV. For 
details of this calculation we refer the reader to Refs. [H, 

23m. 

As mentioned above, a double-counting correction has 
to be applied in order to subtract the contribution to the 
correlation energy already included in the LDA. Several 
forms have been proposed, we will apply the two most 
common approximations, which are the around-mean- 
field (AMF) and the full-localized-limit (FLL) forms, 



E^ MF = U (N - n) - J {N a - n) 

v tr,FLL 
^DC 



U{N- 0.5) — J (N a - 0.5) . 



(6) 
(7) 



where N is the total electronic charge of the impurity 
problem, N a its spin-dependent value, and n the charge 
per spin and orbital. For the corresponding double count- 
ing energies one can find^ 
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with I = 2 the orbital quantum number for 3d electrons. 

Since LSDA calculations give a highly polarized state, 
we perform our spin-polarized DMFT calculations start- 
ing from non-magnetic LDA calculations. 
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FIG. 1: (Color online) The relative total energy of LaFeAsO 
as function of the As height in the unit cell (^-parameter). 
Solid black line (triangles): LDA result. Red line (circles): 
LDA+DMFT using FLL double counting. Blue line (dia- 
monds): LDA+DMFT using AMF double counting. Curves 
are shifted to give similar absolute value of the total energy. 
Vertical dashed line marks the experimental z position. Error 
bars are calculated from averaging several further iterations 
at the self-consistent solution. 



Let us start the discussion of our result with the de- 
termination of the As z position. We did paramag- 
netic LDA+DMFT calculations at inverse temperature 
ft = 40 cV -1 , roughly corresponding to room tempera- 
ture, using the two different types of double counting cor- 
rections mentioned in Sect. [TIJ In Fig. Q] we compare the 
results with the structure optimization within the LDA, 
calculated using the Wien2K package. It is obvious from 
these curves that the inclusion of correlation effects via 
the DMFT significantly improves over the LDA results. 
The ion is pushed away from the iron layers towards 
the experimentally realized z position. Interestingly, the 
choice of the double-counting correction, although having 
almost no effect on the single-particle spectra (see below, 
Sect has some visible effect on the total energy. This 
is most likely due to the very small energy scales that one 
has to deal with in these structure optimizations, where 
already tiny differences can visibly show up. Neverthe- 
less, using the AMF double counting, the As z position 
as determined in LDA calculations (z « 0.634) is cor- 
rected to around z = 0.643, which has to be compared 
with the experimental value of z = 0.651. We attribute 
the larger distance of the As ion from the iron plane to 
the fact that in DMFT calculations the ground state of 
the iron atom is the S=2 high spin state, having slightly 
larger ionic radius then the non-magnetic state realized 
in non-magnetic LDA calculations. 

The small discrepancy that we still see between our 
calculated z values and experimental data is most likely 
due to the neglection of Coulomb interactions between 
the iron and arsenic ions (p-d interactions). It is a very 
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FIG. 2: (Color online) Magnetic properties of LaFeAsO in 
the paramagnetic (high-temperature) phase, calculated using 
AMF double counting corections. Upper panel: Local spin 
susceptibility for inverse temperatures (from top to bottom) 
P — 40, 60, 80, and 100eV _1 , corresponding to temperatures 
T — 290, 193, 145, and 116 K. Lower panel: Static suscepti- 
bility. 



common feature of LDA calculations, that the gap be- 
tween valence bands and ligand bands is too small. Com- 
paring calculated band structure of LaFeAsO with PES 
experiments, the discrepancy is about 1 eV.— In one-shot 
LDA+DMFT calculations this gap can artificially be in- 
fluenced by manually chosen double-counting corrections, 
which is not the case for full self-consistent calculations 
(see below, Sect. The correction of this gap would 
only be possible by the explicit inclusion of p-d interac- 
tions, giving also a repulsion between iron and arsenic 
ions. However, these interactions, without further ap- 
proximations, go well beyond single-site DMFT calcula- 
tions as used here. 

Having established the improved description of the 
crystal structure, we move on to magnetic properties of 
LaFeAsO. From now on, we always use the AMF dou- 
ble counting correction, merit to be more appropriate for 
metallic systems. We use here again the experimental 
value for the As z position in order to compare more 
directly to experiments. In Fig. [2] we plot the local 
susceptibility as function of imaginary time for different 
temperatures. Integrating over imaginary time gives the 
static susceptibility, which is plotted in the lower panel 
of Fig. [2] as black line (open symbols). We see only a 
very weak dependence on temperature, consistent with 
experiments ^ 15 ' 38 Also, the value of x(T) is substan- 
tially enhanced compared to free electrons in agreement 



with Ref. a situation often called 'enhanced' Pauli- 
magnetism. Please note that the upturn at the magnetic 
phase transition in experiments is missed here, since we 
do here calculations only in the paramagnetic phase. 

The instantaneous magnetic moment (equal-time cor- 
relation function) is large, roughly 1.95 /is- How- 
ever, when looking at the ordered moment at low- 
temperatures, the situation is different. To study the 
ordered moment at low temperatures, we performed 
LDA+DMFT calculations allowing for spin polarization 
at T — 116 K, which is well below the magnetic tran- 
sition temperature of TV ~ 140 K. As magnetic order 
pattern we assumed the stripe SDW pattern, as sug- 
gested by LSDA calculations as well as experiments. In 
order to keep calculations as feasible as possible, we used 
ferromagnetic instead of antifcrromagnctic stacking in c- 
direction, but since the distance of the iron layers is very 
large, this approximation is well justified. We used the 
orthorhombic low-temperature unit cell as given in Ref.@, 
with ferromagnetic chains running along the short bonds 
in the xy plane. Doing so, we find an ordered moment of 
m = 0.58 /iB, significantly smaller than our LSDA value 
of m = 1.74 /iB- The value of the magnetic ordered mo- 
ment is almost converged in temperature, since calcula- 
tions for T = 77 K give only slightly larger moments of 
m = 0.60 /iB- In a recent comprehensive LDA+DMFT 
study a value of m = 0.8 /iswas reported for LaFeAsO^i 
the difference in the two results coming from the lar ger 
interaction values U = 5.0 and J = 0.7 used in Ref. [2l| 
(an estimate on the variation of the magnetic moment 
as function of parameters has been given in Ref. Il9r ). 
For comparison, in the first LDA+DMFT study of the 
ordered magnetic moment ^ done for BaFe2As2, a sim- 
ilar reduction of the magnetic moment to m ~ 0.9 /iB 
has been found. The temporal fluctuations, which are 
very strong in the LaFeAsO compound due to its quite 
itinerant nature, hinder the instantaneous moments from 
complete ordering, leaving only a fraction of the moment 
in the ordered state. Similar arguments have been given 
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for the reduction of the moment in Ref. Ojl 

Our findings are in qualitative agreement with a recent 
study on the quenching of the magnetic moment.— How- 
ever, the former study has been done in the paramagnetic 
phase, focusing on the influence of local quantum fluctu- 
ations on the local moment. A direct comparison of the 
values of magnetic moments is therefore not appropriate. 

The reduction of the ordered moment can also be seen 
in the local density of states as shown in Fig. where we 
plot the momentum integrated spectral function for the 
Fe 3d electrons. Real frequency data has been obtained 
by using the stochastic Maximum Entropy method^ 
In LSDA the splitting between majority and minority 
spins is large, whereas we see only a small gap in the 
LDA+DMFT spectra due to the smaller moment. 



IV. CONCLUSIONS 

In this paper, we presented an extension of the previ- 
ously introduced LDA+DMFT approach^ based on the 
augmented linearized plane wave basis to full charge 
self-consistency, and applied this approach to struc- 
tural and magnetic properties of the iron superconduc- 
tor LaFeAsO. We calculated the total energy as function 
of the pnictogen height, and found that the inclusion 
of correlation effects shift the minimum position from 
z = 0.632 to roughly z = 0.644, a much better agree- 
ment with the experimental value of z = 0.651. This 
increased distance of the As ion from the Fe plane is due 
to the high spin state of iron, which is formed due to 
local interactions. 

Considering the magnetic properties, wc calculated the 
local spin susceptibility and found that it shows very 
weak temperature dependence in the paramagnetic state, 
in accordance with enhanced Pauli magnetism. In the 
low-T SDW phase, we calculated the ordered moment 
in the stripe-like antifcrromagnctic phase, and found a 
moment of m » 0.6 /ib, in much better agreement with 
experimental values than the LSDA value, which can be 
(for the experimental crystal structure) as high as 2/is- 

In summary, the inclusion of correlation effects signif- 
icantly improve both structural and magnetic properties 
of LaFeAsO within one set of parameters. This strongly 
points to the importance of local quantum fluctuations 
and correlations for the physics of iron-based supercon- 
ductors. 
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Appendix A: Single-shot versus full charge 
self-consistency 

In this appendix we show that self-consistent calcula- 
tions improve over one-shot calculations regarding the 
choice of the double-counting correction. The most 
straight forward quantity to look at is the impurity self- 
energy on the Matsubara axis, which is not affected by 
any analytic continuation problem to real frequencies. 
Fig.@]shows the result for the iron d xy orbital. It is obvi- 
ous that non-self consistent calculations give a sizable dis- 
crepancy between FLL and AMF double counting, which 
is largely canceled in self-consistent calculations. From 
this plot we can also see, that the FLL one-shot cal- 
culation is in better agreement with the self-consistent 
calculations, whereas AMF is far off. 

Going to the real axis, we can look at the momen- 
tum integrated spectral function, and compare it with its 
LDA result, shown in Fig. [SJ Again, similar as discussed 
above, the agreement between different calculations is 
much better in the self-consistent case (lower panel) , and 
one-shot FLL is again in better agreement. A striking 
difference between non-self-consistent and self-consistent 
calculations is that there is no spurious shift of the As and 
O p states due to the approximate nature of the double 
counting correction. Both calculations show the features 
largely related to As and O at basically the same energy 
as in the LDA calculation. As a result, it is not easily 
possible in self-consistent calculations to use a manually 
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FIG. 5: (Color online) Comparison of the k-summed spec- 
tral function for different LDA+DMFT calculations. Black 
solid lines: LDA result. Green dot-dashed: LDA+DMFT, 
FLL. Red dashed: LDA+DMFT, AMF. Upper panel: non- 
self-consistent. Lower panel: self-consistent calculations. 

adjusted double counting correction for increasing or de- 
creasing the p-d gap. 

Appendix B: LDA+DMFT charge density within 
the (L)APW basis set 

a. Charge density within MT-spheres 

In the LAPW framework the basis functions within 
an MT-sphere a are the radial solutions u"' CT of the 
Schrodingcr equation labeled by the orbital quantum 
number I and spin a and evaluated at a certain lineariza- 
tion energy Ef^ and their corresponding energy deriva- 
tives u[' a evaluated at the same energy. Additional ra- 
dial solutions can be introduced to account for semicore 
states, they are evaluated at a corresponding energy Ef 2 
of the semicore band. The angular and spin dependence 
for the solutions within the MT-spheres is given by cor- 
responding spherical harmonics Y 7 l n (f) and spinors Xa- 
The functions u" ,a {E^ )Yj n Xa, u"'^ (E^Yj^Xa- f° r valence 
and u"' a (E^YjnXa for semicore states will contribute to 
a given eigenvector ^ kv with the corresponding coeffi- 
cients A^(k,a),B^(k, a) and (k, a) , respectively, 
as defined in Ref . y. 

Let us designate the set of these basis functions 
{u^(E^Yixa,u^(E^Yix,,ut' a (E^Yix,} for a 



given MT sphere a and quantum numbers Z, m, a 
as xf'^Yj^Xv an d the set of corresponding coefficients 

i A im ( k > a ) ' ^im (k> CT )' C lm(^ CT )} with which the y con " 
tributes to a given eigenvector VP/^ as Sfi^ (k, a) . Hence, 

within a given MT-sphere the KS eigenvector ^% v (r) 
is expanded as J2i ml S imii k > (X ) x u' T ( r ) Y m(r)x<r, where i 
runs over all radial functions {x} and corresponding co- 
efficients {S}. 

Using those designations, the charge density contribu- 
tion from the states within the energy window W (the 
second term in RHS of Eq. [2]) can be rewritten for a 
given MT-sphere a and for a given spin as 

k IV ij 

X ( Y rn'^Y 
mm' 

x E ^(k,^)(5fe(k,cr))*iV^. (Bl) 

To represent the angular dependence of the charge den- 
sity it is expanded, within a given MT-spheres, in real 
spherical harmonics Y^~(f), 

pZi*)=Y,P™(r)Y™ (?), (B2) 

lm 

For p l Ja(r) = J dfp™ a (r)Y™(f) one obtains from ([Blj) 

k IV ij 

x C lnSimi (B3) 
mm' 

vu'&W 

where C^' imi = I * ^ are the 

corresponding Gaunt coefficients. 



6. Charge density in the interstitial region 

In the interstitial region the LAPW basis functions are 

plain waves — !=e z ( k+G ' )r , where G is the reciprocal lattice 

VV 

vector and V is the unit cell volume, contributing to a 
given KS eigenvector \^% u ) with the corresponding coef- 
ficients Og(k, a). Substituting this into Eq.[5]one obtains 
the contribution of the states within the energy window 
W to the charge density in the interstitial: 
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^a^(k,a)e^ r ^(a^(k,a)) e 

G G' 

£a£(k, CT )e iGr ^(a^(k,cr))*e 



-ik r 



-iG'i 



1 v i/iv' 



1 y vv' 1 



(B4) 



where kc = k + G. 

In the actual computation of (|B4[) one may transform 
the interstitial wave function ^ G a^(k, cr)e iGr to a aux- 
iliary mesh in the real space via the fast Fourier trans- 
form.: 

%(k,a) = 5>£(k,aK GR , (B5) 

G 

therefore getting rid off the double sum over G and G' 
in (|B4[) . The charge density on the axillary mesh then 



reads: 

^( R ) = 7/E E VkMVLMTNt,, (B6) 

k vv'eW 

which is then again transformed back to the reciprocal 
space via inverse FFT: 

/(G)=Eph^ GR - (B7) 

R 
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